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Abstract. We present a straightforward and verified method of deciding whether 
the point x* £ 1", n ^ 1, such that V/(x*) = 0, is the local minimizer, maximizer 
or just a saddle point of a real- valued function /. The method scales linearly with 
dimensionality of the problem and never produces false results. 
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C/^ ! 1. Introduction 



This work is motivated by the practical problem, encountered during 
our studies in physics of magnetic materials. Namely, we wanted to 
investigate properties of simple magnetic systems, consisting of few 
entities, called spins, and treated as classical (i.e. not quantum) 2- or 
J> ' 3-dimensional vectors of unit length. The positions of spins are fixed in 

(N : space (in crystal lattice, for example), but the spins are free to rotate - 

accordingly to the interactions between them and to the strength and 

■ orientation of the external magnetic field. Each geometrical configura- 
tion of the system is characterized by a single number called the free 
energy. It is the so called free energy landscape what we are interested 
in: the positions of (stable) free energy minima, the valleys between 
them and so on. The free energy is a smooth function, defined on the 
open domain spanned by angular variables describing the orientations 
of all the spins involved. The interactions between spins, as well as the 
numbers characterizing the external field (if any), are fixed parameters. 

■ Very similar problem is encountered in computational chemistry, 
where the so called reaction pathways need to be traced. 



2. Standard approach and its deficiencies 

The exploration of the free energy landscape usually begins with solving 
the system of simultaneous equations: 

d n^---,*n) = i = 1>2) ... >n> (1) 

dXj 
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where n is the number of unknowns (variables). 

From now on we will assume that the system (1) has finitely many 
solutions X ^ • X^ j • • • ; X * with p < oo. We will not discuss the potentially 
possible degenerate case, p = oo (countable or not) for purely physical 
reasons: each real system always settles in a state with well defined 
magnetization, at least after sufficiently long time. 

Once the set X* = {x\ : V/(x£) =0, k = 1, 2, . . . ,p} is known, we 
can start the classification procedure. It should tell us which members 
of X* are minimizers, maximizers or correspond to the saddle points of 
/. The usual approach is to investigate the properties of Hessian of the 
function /, calculated for each x\ in turn. Positive definiteness of the 

matrix Hij = 9 ^^dx^dx'^^ x=x l * s * ne sufficient (but not necessary) 
condition for / to have a local minimum at x = x\. Checking whether 
H is positively defined is easy in dimension n = 2 and is routinely 
presented in analytical calculations performed 'by hand', as can be 
seen in many textbooks on magnetism. For n = 2 the process reduces 
to finding whether the expressions 

d 2 f(xj) ^ g/K) and a 2 /(4) 9 2 f«) ( d 2 f(xj) \ 2 

dx\ ' dx\ dx\ dx\ \ dx\dx 2 J 

are all positive (all negative when searching for maximum). 

When the dimensionality of the problem gets higher, the above ap- 
proach becomes more and more tedious, requiring the evaluation of 
many expressions of increasing complexity - determinants of various 
minors of the matrix H (Korn&Korn, 1968). In automated computa- 
tions another approach may appear more efficient, namely finding all 
the eigenvalues of the matrix H (x£). The positiveness (negativeness) 
of all its eigenvalues is also a sufficient condition for H (x* k ) to be 
positively (negatively) defined and, consequently, for x* k to be a local 
minimizer (maximizer) of /. Needless to say that both approaches 
are, except for nearly trivial cases, practically unsuitable for hand 
calculations and we have to rely on computers to perform this task. 

However, both those approaches suffer from two serious problems. 
The first one is inherent to automatic computations, performed with 
limited accuracy. Every investigated point x k , k = 1, . . . ,p, is already 
known only approximately and so is the matrix H (x k ). Rounding errors 
accumulating during either procedure can only worsen this situation 
leading to the unreliable or even false results. 

The second possible deficiency has nothing to do with limited accu- 
racy and is related rather to the properties of the function /. Consider 
for example / (x±, x 2 ) = x\ + x\ having exactly one minimum at x* = 
(0,0). One can easily check, that H (x*) is a singular 2x2 matrix, 
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with the only non- vanishing element d 2 f jdx\ = 2. Since the matrix is 
diagonal then we have immediately its all eigenvalues: Ai = 2, A2 = 
- not all positive. No conclusion concerning x* = (0, 0) is thus possible 
during exact calculations. It is interesting, however, that in automated 
calculations we may arrive at slightly perturbed x* = (0, 5), with <5 7^ 0, 
as a sole candidate for a local minimizer. Now the Hessian is diagonal 
again, with Hu = 2 and H22 = 12x| = 125 2 , leading to different 
conclusions. Depending on the particular value of S, H22 either remains 
equal to zero within the machine accuracy, like before, or is positive. 
For example, working with accuracy of 10 decimal digits we may have: 
5 = 10 -4 and H22 = 2 x 10 -8 > 0, while the relevant component of 
gradient is df/dx2 = 4x| = 4 x 10~ 12 - the number which will be 
rounded down to exactly zero by our computer. Looking at those two 
numbers one is tempted to think that x* = (0, 5) is a true minimizer 
for /. Maybe x* = (0, 0) is another one, for some reason missed by 
gradient-calculating routine. 



3. The interval solution 

Here we present a simple and elegant solution to our problem, based 
on properties of interval calculus. Let us recall the definition of a local 
minimum of a real function / of n variables: 

Definition. We say that / (xi, X2, ■ ■ ■ , x n ) has a local minimum at x* 
when 

3 £ >o Vxev(f) \\x - x*\ < e f (x*) < f (x) , 

where ||-|| is any norm defined in R n , and T>(f) C K" is the domain of 
the function /. 

In simple words: moving away from the point x*, but within the 
limited range e, always leads to the increase of the function value 
compared to f (x*). 

3.1. Wrong, naive test for minimum and why it fails 

Let us try to make direct use of the definition above and let's evalu- 
ate the function / at the following points: (x*, x\, ...,#£ — £,..., x*), 
(x*, X2, ■ ■ ■ , x£ + s, . . . , x*), for some fixed (presumably small) value of 
e and k = 1, 2, . . . , n, in hope of reaching the conclusion concerning the 
character of x* - whether it is a local minimum, maximum or a saddle 
point of /. It is well known that this procedure may only accidentally 
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produce the correct answer. The main reason is that it does not sample 
every possible direction around x*. Last but not least - x* may, and 
usually will, slightly differ from the true location of minimum. 

3.2. Interval test 

The naive test produces incorrect results but, fortunately, we know why. 
Nevertheless its simplicity is so tempting that the idea of improving it 
makes sense. All we have to change is the ability to test the behavior 
of the given function in every possible direction with respect to the 
suspected point. To achieve this goal we need to construct a closed 
surface around x* and simply check the range of / on this surface. 

Here are the necessary steps of the interval-oriented algorithm to 
determine the character of each point x* k € X*, k = 1,2, ...,p, i.e. 
satisfying the equation V/ (x*) = 0: 

1. initialization: set k = 1, 

2. fix the attention at point x\ € X*. Calculate the reference value 
V k = f(x* k ), 

3. determine the distances between x* k and all other members of the 
set X* and discover the shortest one, D k , 

4. set e = D k /2, 

5. generate 2n interval boxes around xjjj with the following properties: 

• the center of each box is an image of the center (midpoint) of 
x k shifted by +e or — e along the consecutive coordinate axes, 

• the size of each box in each direction is 2e, except for the shift 
direction in which the width of box is equal to zero. 

6. evaluate / over each newly created box obtaining the intervals , 
F 1 , F 2 , • • • , F£, F n , 

7. count the events: 

• N : the intervals V k and F* (j = 1, 2, . . . , n, A € {+,-}) 
intersect, 

• _/V>: V). > F* (every real number taken from V k is greater 
than any number form F^), and 

• N < :V k < Ff is true. 

The classification of x\ is following: 
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• iV< = 2n =4* / has a local minimum at x£, 

• iV> = 2n f has a local maximum at x^, 

• iV> • 7V< 7^ =>- there is a saddle point at x* k (inflection point 
ifn=l), 

• otherwise the case is undecided. 

8. set k <— k + 1. If k < p then repeat the procedure, starting from 
step 2 else finish. 



4. Discussion and final remarks 

The sketch of the algorithm makes no clear statement whether the 
elements of the set X* belong to R n or rather to IR n - the set of all n- 
dimensional intervals. For the idea itself, as presented here, it is not an 
issue and both interpretations are almost equally good. This is because 
we don't discuss the ways to obtain the set X*. What we require, how- 
ever, is that X* contains all the solutions of an equation V/ = within 
the domain of interest. This is because we have to be able to precisely 
separate every one of such solution from every other member of X* . 
In machine calculations the really important thing is the knowledge 
of guaranteed bounds for each member of X* and the certainty that 
those solutions are separable, even after their uncertainties are taken 
into account. 

The proposed routine avoids the most important trap of the naive, 
incorrect approach. It effectively samples all the directions around the 
suspected point x* and therefore is in full accord with the definition 
of a local minimum. This is because the trial boxes constructed by the 
algorithm make a complete and 'air-tight' surrounding of the suspected 
point x*. In other words every straight line that passes through x* 
must also necessarily intersect two surrounding boxes. The continuity 
of /, which is differentiable and therefore continuous, assures that our 
algorithm is correct. It is the remarkable property of the interval cal- 
culations: the ability of executing infinite and uncountable number of 
operations in a single step. This feature makes possible to convert the 
naive and essentially wrong algorithm into a powerful and reliable tool. 

It may come as a surprise that our e is rather large, contrary to the 
regular use of this symbol, mostly thought as 'being sufficiently small' 
or 'no matter how small'. We prefer to use e this big for a good reason: 
too small value is dangerous and vulnerable to the other trap, namely 
that x* is inexact. On the other hand the bigger e is the wider can 
be the intervals F* and therefore we may obtain 'undecided' result 
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too often. Our prescription sets the safe upper limit for e rather than 
treats it as the one and only correct value. If e had higher value then 
our surface could contain more than one element of X*. 

In practice the set X* will be determined by interval methods (be- 
cause only those methods guarantee that all the candidates for extrema 
can be found within the domain of interest) and therefore each its 
member will have the form of a small box. Setting e equal to the 
width of such box, or only slightly higher, is the first thing coming 
to the mind. One should not forget, however, that the interval calculus 
usually overestimates the ranges of the functions. For this reason the 
range of / calculated for the single face of such a small box is likely to 
have non-empty intersection with range of / It is even certain, if 
the true minimizer happens to be located at the face of x\ currently 
investigated rather than laying at its midpoint. 

Making the surrounding boxes 'thin' in one direction is a trick to 
circumvent the notorious overestimates of interval enclosures. But it is 
not perfect and larger values of e, somewhere between and half 

of the width of x£ should be used. 'Undecided' members of X* may 

be retried with e' = (e+ iwidth (x%)j /2. Allowing e to be smaller 
that the halved width of x* is dangerous: the uncertainty of x* will 
almost surely produce false results, if ever. The other built-in feature 
of the algorithm is implicit partitioning of the investigated surface into 
2n parts. Doing so we also increase our chances of getting smaller 
overestimations. Of course, in practice the explicit form of / is also 
important - the SUE's (Single Usage Expressions), if possible at all, 
are preferred as usually. The thin trial boxes should be constructed 
with care: their edges have to be rounded outwards. 



Appendix 

It is likely that AP-hardness of many interval algorithms is among 
the key factors preventing their wide dissemination, no matter that 
they also deliver only highest quality, verified results. The algorithm 
presented here is different: its complexity per single candidate scales 
linearly with the dimensionality of the problem. Hence it is able to 
outperform any of its classical counterpart based on matrix opera- 
tions. Moreover, it is simple and makes no implicit use of unfounded 
assumptions, like the one that every minimum can be approximated 
by a quadratic form. In addition, it will never produce false results. 
Our solution is one more example of the old truth: we need algorithms 
designed from the very beginning as interval-oriented. 
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